
#ifndef __RBFFUNCTION_H__
#define __RBFFUNCTION_H__

#include "./mlsdeformer/shape3D.h"
#include "./mlsdeformer/constraints.h"
#include "./mlsdeformer/typedefs.h"

#include <gsl/gsl_blas.h>
#include <gsl/gsl_matrix.h>
#include "solve.h"

   
using namespace std;

//using namespace MLS;

class RBFfunction {
   
public:
   
   RBFfunction();
   ~RBFfunction();
   
   void clear();
   void setting(MLS::Constraints *c);
   
   /// Evaluates the RBF function on a point
   Point3 evalf(const Point3 &p);
   Point3 evalf(const Vector3 &p);
   
   /// Computes the RBF-weights
   void compute_lambdas(MLS::Constraints *c);
   
protected:
   /// Computes the inverse of the RBF-matrix
   void compute_iM();
   
   /// The basis function
   inline Number phi(const Vector3& v);
   
   /// Evaluates a point in a term of the sum
   Number polynomial(const Point3& p, uint i);
   
   /// Define the polynomial degree
   void polynomialDegree(uint i);
   
//    void execute_deformation();
//    Point3 RBFfunction::gradf(const Point3 &p);

public:
   MLS::Constraints     constraints;
   //~ MLS::Shape3D         *model;
   
protected:
   
   gsl_matrix           *iM;
   gsl_matrix           *lambdas;
   gsl_matrix           *RightSide;
   
   unsigned             n;           /// control points number
   unsigned             npolyterms;  /// number of polynomial terms
};

#endif
